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Abstract We present results of our pulsar population 
synthesis of ordinary isolated and millisecond pulsars in 
the Galactic plane. Over the past several years, a pro- 
gram has been developed to simulate pulsar birth, evo- 
lution and emission using Monte Carlo techniques. We 
have added to the program the capability to simulate 
millisecond pulsars, which are old, recycled pulsars with 
extremely short periods. We model the spatial distribu- 
tion of the simulated pulsars by assuming that they start 
with a random kick velocity and then evolve through the 
Galactic potential. We use a polar cap/slot gap model for 
7-ray emission from both millisecond and ordinary pul- 
sars. From our studies of radio pulsars that have clearly 
identifiable core and cone components, in which we fit 
the polarization sweep as well as the pulse profiles in 
order to constrain the viewing geometry, we develop a 
model describing the ratio of radio core-to-cone peak 
fluxes. In this model, short period pulsars are more cone- 
dominated than in our previous studies. We present the 
preliminary results of our recent study and the implica- 
tions for observing these pulsars with GLAST and AG- 
ILE. 
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1 Introduction 

In the very near future, the 7-ray telescopes AGILE and 
GLAST will be launched. It is expected that the number 
of identified 7-ray pulsars will greatly increase. The 7- 
ray pulsar Geminga appears to be radio silent or at least 
very weak. It is not clear whether its radio quiet nature 
is caused by a misalignment of 7-ray and radio beams 
due to viewing geometry, or by intrinsically weak radio 
emission. In the past, polar cap models located the 
7-ray emission along the same open field lines as radio 
emission at similar altitudes. In such a model, the 7-ray 
and radio beams are coaxial with a large overlap, pre- 
dicting a larger number of radio-loud than radio-quiet 
7-ray pulsars. On the other hand, the 7-ray emission re- 
gion in the outer gap models [lOjlSj is located near the 
light cylinder where the 7-ray beam is at a large an- 
gle relative to the radio beam that leads the viewer to 
observe high-energy and radio emission originating from 
opposite poles. The outer gap model predicts a much 
larger number of radio-quiet than radio-loud 7-ray pul- 
sars. Thus the ratio of radio-loud to radio-quiet 7-ray 
pulsars could serve as a discriminating signature to dis- 
tinguish between competing models. 



2 Simulation - Assumptions 

The present simulations are extensions of previous works 
by Gonthier et al. ([7], [8] & [10]) that include both nor- 
mal, isolated and millisecond pulsars from the Galactic 
disk, and the emission of 7-rays within the slot gap/polar 
cap model. Two years ago at the Hong Kong meeting [9] , 
we presented simulations of pulsars from the Gould Belt 
accounting for less than 30% of the EGRET unidentified 
7-ray sources that were previously correlated to the loca- 
tion of the Gould Belt. Recent reassessment of the diffuse 
7-ray background by Casandjian & Grenier [2] suggests 
that the Gould Belt is no longer such a significant source 
of EGRET unidentified 7-ray sources. As a result, we do 
not simulate pulsars from the Gould Belt in this study. 
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For both normal and millisecond (ms) pulsars, we 
assume a birth location in the Galactic disk as given by 
Paczynski [19]. We evolve the neutron stars from their 
birth location to the present in the Galactic potential 
defined by Dehnen & Binney [6]. For normal pulsars, we 
assume a supernova kick velocity distribution of Hobbs 
et al. [13j , a uniform initial period distribution from to 
500 ms, initial magnetic field distributions with a decay 
constant as given in Gonthier et al. [8J and a uniform 
birth rate back to 1 Gyr. 

For ms pulsars, we begin the evolution by using the 
magnetic field and supernova kick velocity distributions 
of Cordes & Chernoff [4j . Recent studies of low-mass X- 
ray binary systems (LMXBs) have been able to deter- 
mine the spins of the accreting neutron stars allowing 
for an estimate of their magnetic fields. Lamb & Yu [T7] 
conclude that the properties of these LMXBs are consis- 
tent, if they have magnetic fields between 3 x 10^ G and 
3 X 10® G and accretion rates ranging from the Eddington 
critical rate Me to 3 x lO^^M^. These different accretion 
rates result in different birth lines in the P — P diagram 
(see figure 4 in Lamb & Yu ^7\). We have included an 
approximate procedure to take into account this distri- 
bution of accretion rates by dithering the intercept of the 
birth line described as 

log(P) = ^ log(F) - 14.9. (1) 

While we tried a Gaussian distribution of the birth 
line, we obtained better agreement with a uniform distri- 
bution of birth lines between Me to 3 x IO^^MbUI]- We 
assume a uniform birth rate for ms pulsars back to 12 
Gyr. We explore various power laws of the magnetic field 
distribution getting better agreement with n{B) oc B^^ 
and with a Bmm = 2 x lO^G. Given the magnetic field, 
a selected birth line and the pulsar age, we obtain the 
present period and period derivative. 

We use a larger Galactic scale height of 200 pc for 
ms pulsars, compared to the one used for normal pulsars 
of 75 pc flW . Since the supernova kick velocities of ms 
pulsars are much smaller than those of normal pulsars, 
most ms pulsars remain bound to the Galaxy, oscillating 
in and out of the Galactic plane with time. We evolve a 
large number of neutron stars generated in our Monte 
Carlo simulation to determine the equilibrium spatial 
distribution of ms pulsars, which we then use in subse- 
quent simulations. A scale height of 410 pc of the evolved 
ms pulsars is in good agreement with a scale height of 
500 pc (exponential scale) of Cordes & Chernoff ^4^ and 
with a scale height of 410 pc of LMXBs pi] . 

For the radio luminosity, we assume that radio pul- 
sars are standard candles and follow the prescription of 
Arzoumanian, Chernoff & Cordes [1] (ACC) given by the 
expression 

L = 2.1 X 10^^P"^-^P° '*mJy • kpc^ • MHz. (2) 

However, to obtain good agreement between the mea- 
sured and simulated flux and distance distributions, we 



needed to reduce the luminosity by a factor of 73 and 
200 for normal pulsars and ms pulsars, respectively. We 
follow the geometric model of the core beam in ACC [1] 
and the cone model of Kijak & Gil [TH] described in this 
meeting by Harding, Grenier, & Gonthier [12]. In the 
study of ACC, the ratio of the core-to-cone peak fluxes 
had a dependence, resulting in a dominance of the 
core component for short period pulsars. As discussed 
by Harding et al. [12], we find compelling evidence that 
short period pulsars are more cone dominated than in 
the ACC model and have adopted a different model to 
characterize this relationship (see [12]). 

The "detection" of our simulated radio pulsars is ac- 
complished using the characteristics of ten radio surveys, 
six of which are at a frequency around 400 MHz - Arecibo 
3, Arecibo 2, Greenbank 2, Greenbank 3, Molongo 2 
and Parkes 2, and four surveys are around 1400 MHz 
- Parkes 1, Jodrell Bank 2, Parkes Multibeam, and we 
have recently added the Swinburne Intermediate Lati- 
tude survey at 1400 MHz. These surveys provide us with 
1208 normal and 50 ms pulsars to which we normal- 
ize our simulation. We adjust the radio luminosity us- 
ing only the Parkes MB pulsar survey in order to get 
a neutron star birth rate of about 2 per century. This 
survey has over 800 detected pulsars, and we have the 
best description of the minimum radio flux 5*111111 of the 
ten surveys. Normalization to 1208 radio pulsars seen 
by the group of 10 surveys is then the only overall ad- 
justment made, thereby, allowing the prediction of birth 
rates, radio-loud and radio-quiet 7-ray pulsars detected 
by the instruments EGRET, AGILE and GLAST. Simu- 
lated neutron stars whose radio flux is below the survey 
threshold, S'min, are assumed to be radio-quiet. 

The simulation of the 7-ray emission occurs over two 
regions in the P — P diagram separated by the curva- 
ture radiation pair death line (CRPDL), below which the 
curvature radiation 7 rays no longer produce electron- 
positron pairs. Above this death line, emission originates 
from low altitude pair cascades on the inner edge of the 
slot gap ^8j| , as well as from primaries accelerating in the 
slot gap at high altitude forming a caustic component 
[T^ . Below this death line, the slot gap dissolves due to 
the less effective screening of the electric fleld, leading to 
extended emission over the entire polar cap. The 7-ray 
emission above CRPDL is most important for normal 
pulsars while emission below the CRPDL is the domi- 
nant 7-ray emission mechanism for ms pulsars. 

The simulated 7-ray flux is compared to the all sky 
threshold maps for EGRET, AGILE and GLAST. We 
include the revised EGRET map that includes the dark 
clouds (see [2]) and the revised GLAST map (after DC2) 
without the dark clouds (Grenier private communica- 
tion). The all sky map for AGILE (Pellizzoni private 
communication) has not been recently updated. 
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Fig. 1 Distributions of various characteristics of normal pulsars indicated as detected pulsars (shaded histograms) and 
simulated (unshaded histograms) pulsars from the Galactic plane. 



3 Results 

In Figure [U we compare simulated distributions (un- 
shaded histograms) of various characteristics of normal 
pulsars with those detected (shaded histograms) by the 
select group of 10 radio surveys. For the flux distributions 
at 400 and 1400 MHz, we have used the values given in 
the ATNF pulsar catalog for the detected distributions. 
In the simulation, we assign a flux at 400 MHz if the 
pulsar is detected by one of the low frequency surveys 
in our group and likewise for the fluxes at 1400 MHz. 
As can be seen in the comparisons of these distributions, 
the model simulation over-predicts the number of young 
pulsars with short periods. The distance distribution is 
now much improved from our previous results ^10. mainly 
due to the inclusion of the actual width of the simulated 
pulse profile of each pulsar in the calculation of the sur- 
vey flux threshold, 5min- Overall we generally see good 
agreement between the simulation and detected distri- 
butions with a predicted birth rate of 1.8 normal pulsars 
per century. 

The simulation again overestimates the number of 
young and distant ms pulsars as shown in Figure 2. Evi- 
dent in the histograms is the limitation of a sample of 50 
ms pulsars detected by the select group of ten radio sur- 



veys. We simulate 500 ms pulsars and then normalize to 
the detected 50 in order to obtain smoother histograms 
for the simulated pulsars. The simulation again overes- 
timates young ms pulsars that are more distance than 
those detected. The 400 MHz surveys in the simulation 
seem to be a bit too sensitive. The predicted birth rate 
is 6.8 X 10""* per century, which is very close to the birth 
rates estimated by Lorimer [18] of 2.9 x 10~^ per century 
and by Kiel & Hurley [TS] of 6.5 x 10^* per century of 
neutron stars from LMXBs. 

In Figure 3, we present the Aitoff plots and P — P di- 
agrams for detected (left) and simulated (right) normal 
(dots) and MS (crosses) pulsars. Both simulated normal 
and ms pulsar distributions in the Galaxy are similar to 
the distributions of those detected. As ms pulsars are 
closer than normal pulsars, they appear with a larger 
out-of-plane distribution. We are not quite able to repro- 
duce the "bunching" of ms pulsars with periods between 
2 and 10 ms. However, the uniform distribution of the 
birth line as suggested by Lamb & Yu [T7j does seem to 
be important in order to get fairly good agreement in the 
P — P diagram with detected ms pulsars. The distribu- 
tion of detected normal pulsars in the P — P diagram is 
more or less reproduced by the simulation, but the sim- 
ulated distribution is too broad in period with too few 
high- field pulsars. We do find that magnetic field decay 
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Fig. 2 Distributions of various cliaracteristics of millisecond pulsars indicated as detected (shaded histograms) and simulated 
(open histograms) pulsars from the Galactic plane. 



with a decay constant of 2.8 Myr is necessary for normal 
pulsars in order to reproduce the detected distribution. 
We did not incorporate field decay in the simulation of 
ms pulsars. 

In Table 1, we present the simulated normal and ms 
7-ray pulsar statistics for radio-loud and radio-quiet 7- 
ray pulsars for various instruments indicated as well as 
the ones detected by EGRET. For normal pulsars, we 
show separately the number of 7-ray pulsars having de- 
tected emission from the low altitude and high altitude 
slot gap, as both of these radiation mechanisms con- 
tribute above the curvature radiation pair death line. As 
expected the low altitude slot gap beam is more co-axial 
with the radio beam geometry and results in a greater 
number of radio-loud 7-ray pulsars. On the other hand, 
the high altitude emission of the slot gap occurs along 
the last open field lines all the way out to the light cylin- 
der, at which point the 7-ray beam makes a larger angle 
to the radio beam and more radio-quiet 7-ray pulsars 
are detected. In this respect, the high altitude emission 
of the slot gap resembles the correlations of radio and 
7-ray beams of outer gap emission. 

The ms 7-ray pulsars in the P — P diagram appear 
below the curvature radiation pair death line and, there- 
fore, their 7-ray emission arises from the emission over 



the entire polar cap region as the slot gap geometry dis- 
appears near the curvature death line. 



4 Conclusions 

We present the preliminary results of a population syn- 
thesis study of both normal and millisecond pulsars. We 
include a radio beam geometry that includes a single core 
and a single cone beam with a new dependence on the 
core-to-cone peak fiuxes. In this new model, as discussed 
in this conference by Harding et al. [12], short period 
pulsars are more cone dominated than in our previous 
studies that follow the ACC model [1]. This relation- 
ship is especially important for millisecond pulsars. We 
use the same radio model for both normal and millisec- 
ond pulsars. We describe the radio luminosity using the 
prescription of ACC (Equation 2), but with luminosity 
decreased by a factor of 73 and 200 for normal and ms 
pulsars, respectively, to achieve reasonable neutron star 
birth rates, flux and distance distributions. The same set 
of parameters are then used in the simulation of normal 
and millisecond pulsars. The simulations are run until 
the same number of normal or millisecond radio pulsars 
are simulated as detected by a select group of ten radio 
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Fig. 3 Aitoff plots of normal (dots) and millisecond (crosses) pulsars detected (left) by the select group of ten radio surveys 
and simulated (right). Lower plots are the P — P diagrams of normal (dots) and millisecond (crosses) pulsars detected (left) 
and simulated (right). Dashed lines represent the pair death lines for curvature radiation (CRPDL) and for nonresonant 
inverse Compton scattering (NRICSPDL). Solid lines represent the traditional magnetic surface field strength, assuming a 
constant dipole spin-down field. 

Table 1 Simulated radio-loud and radio-quiet 7-ray pulsar statistics 

Normal Pulsars Millisecond Pulsars 
low + high altitude 

Instrument Radio- Loud Radio-Quiet Radio-Loud Radio- Quiet 

EGRET Detected 6 1 1 

EGRET Simulated 25+2 10+22 1 2 

AGILE Simulated 51+1 16+53 2 3 

GLAST Simulated 94+6 51+100 4 10 



surveys providing an overall normalization of the sim- 
ulation. With the same set of parameters for the radio 
beam geometry and luminosity, we achieve very reason- 
able agreement with the detected distributions of various 
pulsar characteristics for both normal and millisecond 
pulsars. This conclusion is quite remarkable. 

However, we expect the radio beam geometry of mil- 
lisecond pulsars to be significantly different, especially 
for pulsars with periods less than 100 ms, due to the spe- 
cial relativistic effects of aberration, time delays and the 
sweepback of the magnetic field. These effects and their 



contributions are discussed in this meeting by Harding 
et al. [12]. We hope to soon include this component to 
the radio emission in our Monte Carlo code to more ad- 
equately describe millisecond pulsars. 

In previous studies, we only included low altitude 7- 
ray emission from the slot gap/polar cap model. Both 
curvature and synchrotron radiation contribute along the 
last open field line about three stellar radii above the 
surface resulting in a conical beam symmetric about the 
magnetic axis. As a result, the 7-ray beam is strongly 
correlated to the core and cone radio beams assumed 



6 



Peter L. Gonthier et al. 



in our simulations. Such a correlation is reflected in the 
larger number of radio-loud than radio-quiet 7-ray pul- 
sars with a ratio of ratio-loud to radio-quiet of 1.8. As 
discussed by Harding et al. [T^] at this meeting, the elec- 
tric field parallel to the magnetic field along the very last 
open field lines is not screened, forming a especially nar- 
row slot gap in the case of short period pulsars. The 7-ray 
beams at high altitude are concentrated (caustic) along 
the last open field lines with the beam being at large 
angle to the radio emission. High altitude 7-ray emission 
then leads to a larger number of radio-quiet 7-ray pulsars 
being simulated by this mechanism with a ratio of radio- 
loud to radio-quiet of about 0.06 for this component as 
indicated in Table 1. The viewing geometry defining the 
impact angle will determine the respective contributions 
of the low and high altitude emission. These components 
will have different signatures in their pulse profiles. 

EGRET saw pulsed emission from only one 7-ray 
ms pulsar, J0218-I-4342. Though not very optimistic, the 
simulations predict that GLAST should detect of the 
order of 10 ms pulsars as point sources. Of the 10 radio- 
quiet 7-ray ms pulsars only 1 is expected to be detected 
through blind searches. 

Our study of normal, isolated pulsars indicates that 
GLAST will be able to detect about 100 radio-loud and 
151 radio-quiet 7-ray pulsars. These numbers are sig- 
nificantly lower than presented in Gonthier et al.[8], 
primarily because the new GLAST point-source detec- 
tion threshold is higher after the great detailed study 
of the GLAST second data challenge (DC2). The chal- 
lenge for GLAST will be to detect the radio-quiet 7- 
ray pulsars through its ability to perform blind searches 
with those pulsars whose 7-ray fluxes are higher than 
lO"'^ photons/ (cm^ • s) (Grenier private communica- 
tion). Out of the 151 radio-quiet 7-ray pulsars GLAST 
should be able to detect pulsations in 50 pulsars. The 
simulation modeling the high and low altitude 7-ray 
emission within the theoretical framework of the polar 
cap model suggests that the expected ratio of radio-loud 
to radio-quiet 7-ray pulsars is about 1.5, while the outer 
gap models predict many more radio-quiet 7-ray pulsars 
and a ratio of about 0.1 radio- loud to radio-quiet 7-ray 
pulsars is expected JT]. Note that the outer gap ratio is 
similar to the ratio for the high altitude emission alone. 
Further detailed studies of the correlations between the 
light curves of radio and 7-ray pulses may provide the 
best signature to distinguish between the outer gap and 
polar cap models. 
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